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1.  EXECUTIVE  SUMMARY. 


1.1.  The  Planning  Model  based  on  Projection  Methodology  (PM2). 

The  U.S.  Army  Materiel  Systems  Analysis  Activity  (AMSAA)  has  developed  a 
new  reliability  growth  planning  model.  The  new  model  is  referred  to  as  the  Planning 
Model  based  on  Projection  Methodology,  or  PM2.  PM2  can:  (1)  aid  in  constructing  a 
reliability  growth  planning  curve  over  a  developmental  test  program  useful  to  program 
management;  (2)  serve  as  a  baseline  against  which  reliability  assessments  can  be 
compared  and;  (3)  highlight  the  need  to  management  when  reallocation  of  resources  is 
necessary. 

1.2.  PM2  Assumptions. 

There  are  a  number  of  reasonable  assumptions  associated  with  PM2.  The  first 
assumption  is  that  there  are  a  large  number  of  potential  failure  modes.  As  a  rule  of 
thumb,  the  potential  number  of  failure  modes  should  be  at  least  five  times  the  number  of 
failure  modes  that  are  expected  to  be  surfaced  during  the  planned  test  period.  Second, 
each  failure  mode  time  to  first  occurrence  is  assumed  exponential.  Finally,  it  is  assumed 
that  each  failure  mode  occurs  independently  and  causes  system  failure. 

1.3.  PM2  Limitations. 

All  reliability  growth  planning  models  have  limitations.  The  first  limitation 
associated  with  PM2  is  that  the  portion  of  testing  utilized  for  reliability  growth  planning 
should  be  reflective  of  the  Operation  Mode  Summary/Mission  Profile  (OMS/MP).  This 
limitation  is  not  unique  to  PM2  -  it  is  a  limitation  associated  with  all  reliability  growth 
planning  models.  Second,  one  must  postulate  a  baseline  test  schedule.  That  is,  one  must 
determine  the  number  of  hours,  or  miles,  per  unit  per  month  over  the  test  period.  Finally, 
fix  implementation  periods  must  be  specified  within  the  planned  test  schedule. 

1.4.  PM2  Benefits. 

There  are  a  number  of  benefits  associated  with  the  new  planning  model.  PM2  is 
unique  in  comparison  to  other  reliability  growth  planning  models  in  that  it  utilizes 
planning  parameters  that  are  directly  influenced  by  program  management.  Some  of  these 
parameters  include:  (1)  the  initial  system  MTBF;  (2)  the  fraction  of  the  initial  failure  rate 
addressable  via  corrective  action  (referred  to  as  management  strategy);  (3)  the  goal 
system  MTBF;  (4)  the  average  fix  effectiveness  of  corrective  actions;  (5)  the  duration  of 
developmental  testing  and;  (6)  the  average  delay  associated  with  fix  implementation.  A 
second  benefit  of  PM2  is  that  the  model  can  detennine  the  impact  of  changes  to  the 
planned  test  schedule,  and  associated  fix  implementation  periods.  Third,  PM2’s 
measures  of  programmatic  risk  are  not  sensitive  to  the  length  of  the  initial  test  phase 
(which  is  a  limitation  of  the  MIL-HDBK-189  planning  model).  Finally,  PM2  can  be 
applied  to  programs  with  limited  opportunities  for  implementation  of  corrective  actions. 
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1.5.  Overview. 


In  following  sections  of  this  report,  exact  expressions  for  the  expected  number  of 
surfaced  failure  modes  and  system  failure  intensity  as  functions  of  test  time  are  presented 
under  the  assumption  that  the  surfaced  modes  are  mitigated  through  corrective  actions. 
These  exact  expressions  depend  on  a  large  number  of  parameters.  Functional  forms  are 
derived  to  approximate  these  quantities  that  depend  on  only  a  few  parameters.  Such 
parsimonious  approximations  are  suitable  for  developing  reliability  growth  plans  and 
portraying  the  associated  planned  growth  path.  Simulation  results  indicate  that  the 
functional  fonn  of  the  derived  parsimonious  approximations  can  adequately  represent  the 
expected  reliability  growth  associated  with  a  variety  of  patterns  for  the  failure  mode 
initial  rates  of  occurrence.  A  sequence  of  increasing  MTBF  target  values  can  be 
constructed  from  the  parsimonious  MTBF  projection  approximation  based  on:  (1) 
planning  parameters  that  determine  the  parsimonious  approximation;  (2)  corrective  action 
mean  lag  time  with  respect  to  implementation  and;  (3)  the  test  schedule  that  gives  the 
number  of  planned  Reliability,  Availability,  and  Maintainability  (RAM)  test  hours  per 
month  and  specifies  corrective  action  implementation  periods. 
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2.  INTRODUCTION. 


2.1.  Background. 

To  mature  the  reliability  of  a  complex  system  under  development  it  is  important 
to  formulate  a  detailed  reliability  growth  plan.  One  aspect  of  this  plan  is  a  depiction  of 
how  the  system’s  reliability  is  expected  to  increase  over  the  developmental  test  period. 
The  depicted  growth  path  serves  as  a  baseline  against  which  reliability  assessments  can 
be  compared.  Such  baseline  planning  curves  for  Department  of  Defense  (DoD)  systems 
have  frequently  been  developed  in  the  past  utilizing  the  assumed  reliability  growth 
pattern  specified  in  Military  Handbook  189  (MIL-HDBK-189)  (Department  of  Defense, 
1981).  This  growth  relationship  is  between  the  reliability,  expressed  as  the  mean  test 
duration  between  system  failures  and  a  continuous  measure  of  test  duration  such  as  time 
or  mileage.  The  equation  governing  this  growth  pattern  was  motivated  by  the  empirically 
derived  linear  relationship  observed  for  a  number  of  data  sets  by  Duane  (1964),  between 
the  developmental  system  cumulative  failure  rate  and  the  cumulative  test  time  when 
plotted  on  a  log-log  scale. 

2.2.  Purpose. 

In  this  paper  we  obtain  a  non-empirical  relationship  between  the  mean  test 
duration  between  system  failures  and  cumulative  test  duration  that  can  be  utilized  for 
reliability  growth  planning.  This  relationship  is  derived  from  a  fundamental  relationship 
between  the  expected  number  of  failure  modes  surfaced  and  the  cumulative  test  duration. 
For  convenience,  we  shall  refer  to  the  test  duration  as  test  time  and  measure  the  reliability 
as  the  mean  time  between  system  failures  (MTBF).  The  functional  form  of  this 
fundamental  relationship  is  well  known  and  is  easily  established  without  recourse  to 
empiricism  (Crow,  1982).  We  obtain  an  approximation  to  this  relationship  that  is 
suitable  for  reliability  growth  planning.  One  significant  advantage  to  our  approach  is  that 
it  does  not  rely  on  an  empirically  derived  relationship  such  as  the  Duane  based  approach. 
We  shall  show  how  the  cumulative  relationship  between  the  expected  number  of 
discovered  failure  modes  and  the  test  time  naturally  gives  rise  to  a  reliability  growth 
relationship  between  the  expected  system  failure  intensity  and  the  cumulative  test  time. 
The  presented  approximation  for  the  resulting  growth  pattern  avoids  a  number  of 
deficiencies  associated  with  the  Duane/MIL-HDBK-189  approach  to  reliability  growth 
planning. 

2.3.  Study  Overview. 

In  Section  3  we  highlight  a  number  of  issues  associated  with  the  Duane/MIL- 
HDBK-189  approach  to  reliability  growth  planning.  Section  4  develops  the  exact 
expected  system  failure  intensity  and  parsimonious  approximations  suitable  for  reliability 
growth  planning.  These  functions  of  test  time  are  derived  from  the  exact  and  planning 
approximation  relationships  between  the  expected  number  of  surfaced  failure  modes  and 
the  cumulative  test  time.  The  exact  relationship  is  expressed  in  terms  of  the  number  of 
potential  failure  modes,  k,  and  the  individual  initial  failure  mode  rates  of  occurrence. 
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Parsimonious  approximations  to  this  relationship  are  obtained.  The  first  approximation 
utilizes  k  and  several  additional  parameters.  The  second  approximation  discussed  is  the 
limiting  form  of  the  first  approximation  as  k  increases.  This  approximation  is  suitable  for 
complex  systems  or  subsystems.  The  approximations  are  derived  through  consideration 
of  an  MTBF  projection  equation.  This  equation  arises  from  considering  the  problem  of 
estimating  the  system  MTBF  at  the  start  of  a  new  test  phase  after  implementing 
corrective  actions  to  failure  modes  surfaced  in  a  preceding  test  phase.  This  MTBF 
projection  has  been  documented  in  (Ellner  et  ah,  2004)  and  is  described  in  Section  4. 

Section  5  contains  simulation  results.  The  simulations  are  conducted  to  obtain 
actual  patterns  for  the  cumulative  number  of  surfaced  failure  modes  versus  test  time  for 
random  draws  of  initial  mode  failure  rates  from  several  parent  populations,  and  for  a 
geometric  sequence  of  initial  mode  failure  rates.  The  resulting  stochastic  realizations  are 
compared  to  the  theoretical  expected  number  of  potential  surfaced  failures  modes  and  to 
the  parsimonious  approximations.  Random  draws  for  mode  fix  effectiveness  factors 
(FEFs)  (fraction  reductions  in  initial  failure  mode  rates  of  occurrence  due  to  mitigation) 
are  used  to  simulate  corrective  actions  to  surfaced  failure  modes.  Using  the  simulated 
corrective  actions,  the  relationship  between  the  expected  system  failure  intensity  and 
cumulative  test  time  is  simulated  for  various  sets  of  mode  initial  failure  rates.  This 
relationship  is  obtained  under  the  assumption  that  the  system  failure  intensity  associated 
with  a  cumulative  test  time  t  reflects  implementation  of  corrective  actions  to  the  modes 
surfaced  by  t  with  the  associated  randomly  drawn  FEFs.  The  resulting  system  MTBF 
versus  test  time  relationship  is  compared  to  the  corresponding  relationship  established  for 
planning  purposes. 

Section  6  derives  expressions  for  a  reliability  projection  scale  parameter  that  is 
utilized  in  the  parsimonious  approximations.  The  projection  parameter  is  expressed  in 
tenns  of  basic  planning  parameters.  The  resulting  MTBF  approximations  are  compared 
to  the  reciprocals  of  the  exact  expected  system  failure  intensity  and  stochastic  realizations 
of  the  system  failure  intensity,  and  to  MIL-HDBK-189  MTBF  approximations  based  on 
planning  parameters.  The  comparisons  are  done  for  several  reliability  growth  patterns. 

Section  7  addresses  the  relationship  between  the  theoretical  upper  bound  on  the 
achievable  system  MTBF,  termed  the  growth  potential,  and  the  planning  parameters.  The 
projection  scale  parameter  considered  in  Section  6  is  then  expressed  in  tenns  of  planning 
parameters  and  the  MTBF  growth  potential.  It  is  shown  that  the  scale  parameter  becomes 
unrealistically  large  if  the  goal  MTBF  is  chosen  too  close  to  the  growth  potential  or  if  the 
allocated  test  time  to  grow  from  the  initial  to  goal  MTBF  is  inadequate. 

Section  8  indicates  how  to  construct  a  sequence  of  MTBF  target  values  that  start 
at  an  expected  or  measured  initial  MTBF  and  end  at  the  goal  MTBF.  It  is  shown  that  the 
parsimonious  approximation  to  the  reciprocal  of  the  expected  system  failure  intensity  can 
be  used  for  this  purpose  in  conjunction  with  a  test  schedule  that  specifies  the  expected 
monthly  RAM  hours  to  be  accumulated  on  the  units  under  test  and  the  planned  corrective 
action  periods. 
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3.  MIL-HDBK- 189  IDEALIZED  GROWTH  CURVE. 


3.1.  Background. 

The  frequently  referenced  United  States  Department  of  Defense  MIL-HDBK- 189 
(Department  of  Defense,  1981)  utilizes  an  idealized  reliability  growth  pattern  with  failure 
intensity  function  pM  (t)  given  by, 

pM{i)  =  XMPMt^  (1) 

Letting  N M  (0  denote  the  number  of  failures  by  t,  the  corresponding  expected  number  of 
failures  by  t  is  given  by, 

E{NM(t))  =  XMt^  (2) 

For  planning  purposes,  the  handbook  divides  the  allocated  test  time  T  into  p  test 
phases  that  end  at  the  cumulative  test  times  q  <t2  < ...  <  tp  =  T .  Reliability  growth  may 

occur  within  a  test  phase.  However,  the  MIL-HDBK- 189  approach  does  not  assume  the 
growth  pattern  governed  by  Equation  (1)  holds  within  each  test  phase.  For  example,  no 
corrective  actions  may  be  applied  within  some  test  phases.  The  test  phases  typically  are 
separated  by  blocks  of  calendar  time  during  which  a  significant  number  of  corrective 
actions  are  implemented  to  failure  modes  surfaced  in  the  preceding  test  phase.  Thus 
jumps  in  reliability  are  typically  expected  from  test  phase  to  test  phase.  For  planning,  the 
handbook  only  assumes  that  Equation  (2)  holds  at  the  ends  of  the  test  phases,  i.e.,  for 

t  =  tf.  In  particular,  the  points  with  coordinates  ti ,  are  assumed  to  lie  on  a 

V  6 

straight  line  on  a  log-log  plot  of  versus  t,  with  slope  equal  to  -aM .  The  value 

h 

aM  is  called  the  growth  rate.  The  parameter  pM  that  appears  in  Equation  (2)  is  equal  to 
1  -aM.  The  assumption  in  the  handbook  is  motivated  by  Duane’s  empirical  relationship 
(Duane,  1964).  Duane  observed  that,  for  a  number  of  data  sets,  the  logarithm  of  the 
cumulative  failure  rate  versus  the  logarithm  of  the  cumulative  test  time  tended  to  display 
a  linear  relationship.  Note,  however,  Duane’s  observations  were  based  on  fitting  straight 
lines  to  test  data  using  a  log-log  scale  and  thus  presumably  applied  to  the  observed 
pattern  of  growth  during  a  test  phase  as  opposed  to  across  test  phases  although  this  is  not 
addressed  in  (Duane,  1964).  In  fact,  although  for  planning  purposes  MIL-HDBK- 189 
does  not  assume  Equation  (2)  holds  within  all  the  test  phases,  the  handbook  does  utilize  a 
power  law  relationship  for  estimating  reliability  in  a  test  phase,  i.e.,  for  tracking 
reliability  growth  within  a  test  phase.  The  use  of  Equation  (2)  for  this  purpose  seems 
more  tied  to  Duane’s  observations. 

3.2.  Depicting  the  Planned  Reliability  Growth. 

The  MIL-HDBK- 189  approach  utilizes  average  MTBF  values  over  each  test 
phase  to  depict  the  planned  reliability  growth  path.  For  test  phases  during  which  no 
corrective  actions  are  expected  to  be  implemented  this  average  MTBF  would  also  be  the 
instantaneous  MTBF.  The  one  structured  method  presented  in  the  handbook  to  obtain 
these  test  phase  average  MTBF  values  is  based  on  first  specifying  an  idealized  curve  that 
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satisfies  Equation  (2)  at  the  end  of  each  test  phase.  The  idealized  curve  is  frequently 
determined  by  specifying  a  planned  or  assessed  average  MTBF,  Mx ,  over  the  initial  test 
period,  the  total  test  time  over  all  the  test  phases,  T,  and  the  goal  MTBF,  Mc ,  to  be 
attained  at  T.  The  growth  rate  aM  can  then  be  solved  for  and  the  constants  pM  =  1  -  aM 

and  /.\,i  that  appear  in  Equation  (2)  can  be  obtained.  For  test  phase  i,  the  average  failure 
rate  <p;  is  defined  by, 

^  _  E{Nu(tt))-E{NMqtj)  (3) 

h  -  h-\ 

The  corresponding  average  MTBF  for  phase  i  is  taken  to  be  the  reciprocal  of  p .  The 
pattern  of  test  phase  MTBF  averages  so  obtained  do  not  explicitly  take  into  account 
parameters  that  can  be  directly  influenced  by  program  management.  Such  parameters 
include  the  management  strategy  (MS),  the  average  fix  effectiveness  factor  ( /irl ),  the 
corrective  action  lag  time,  and  the  scheduled  monthly  RAM  test  hours  and  corrective 
action  periods.  If  one  explicitly  takes  into  account  these  factors,  the  growth  pattern 
exhibited  by  the  resulting  test  phase  MTBFs,  even  when  displayed  on  a  test  time  basis, 
will  often  look  more  irregular  than  the  pattern  of  average  test  phase  MTBFs  obtained 
from  the  MIL-HDBK-189  approach.  In  particular,  this  growth  pattern  may  not  be  well 
represented  by  a  pattern  consistent  with  Equation  (1). 

The  reciprocal  of  the  idealized  failure  intensity  given  in  (1)  is  considered  to  be  a 
representation  of  the  overall  MTBF  growth  trend  over  the  test  program  after  the  first  test 
phase.  Note  for  growth  one  has  aM  greater  than  zero.  Thus  as  t  approaches  zero  the 
MTBF  implied  by  Equation  (1)  approaches  zero.  Therefore,  for  planning  purposes,  the 
handbook  utilizes  Equation  (1)  to  represent  the  overall  growth  trend  only  for  t>tl  .  The 


handbook  simply  utilizes  a  constant  or  average  failure  rate,  </\  =  M  ,  over  the  first  test 
phase.  The  constant  <j>x  is  chosen  such  that  Equation  (2)  is  satisfied  for  t  =  tx.  Doing  so,  it 
follows  that  the  MTBF  growth  trend  consistent  with  Equation  (1)  for  t>tx  and  </>x  is  given 


by, 


MTBF{t)= 


/  V«A/ 


0  <t<tx 
t>tx 


(4) 


3.3.  MIL-HDBK-189  Planning  Model  Issues. 


In  using  Equation  (4)  one  must  be  careful  not  to  automatically  equate  Mx  to  the 
planning  parameter  M,,  defined  as  the  initial  MTBF.  In  general,  M,<MX.  The  two 
MTBFs  should  be  equated  only  if  no  growth  is  planned  over  the  first  test  phase,  since  Mx 
is  the  planned  average  MTBF  over  the  initial  test  phase. 

The  growth  rate  aM  is  used  as  a  measure  of  programmatic  risk  with  respect  to  being 
able  to  grow  from  Mx  to  Ma  =  MTBF(T)  in  test  time  T.  The  higher  the  aM  relative  to  past 
experience  the  greater  the  risk  of  attaining  Ma.  From  Equation  (4)  we  can  see  that 
MTBF(t  )  is  a  strictly  increasing  function  of  the  ratio  T/tx  and  can  be  made  as  large  as 
desired  by  making  tx  sufficiently  small.  Thus  for  any  given  T,  M  ,  and  growth  rate  aM 
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one  can  always  find  a  small  enough  /,  such  that  MTBF(T )  will  equal  the  desired  value. 
This  implies  that  aM  as  a  measure  of  programmatic  risk  is  only  as  meaningful  as  the 
choice  of  tx .  In  particular,  one  should  guard  against  artificially  lowering  aM  by  selecting 
tx  so  small  that  no  significant  amount  of  fix  implementation  is  expected  to  occur  until 
during  a  corrective  action  period  that  is  beyond  /, .  The  strong  dependence  of  the  global 
parameter  aM  on  the  length  of  the  initial  test  phase  is  not  a  desirable  attribute  for 
planning  purposes. 

Finally,  we  note  that  Equation  (4)  implies  that,  even  with  a  reasonable  choice  for 
tx ,  any  value  of  Mc  can  eventually  be  obtained  since  there  is  no  upper  limit  implied  by 
Equation  (4).  This  is  true  even  using  a  growth  rate  that  appears  to  be  reasonable  based  on 
past  experience  with  similar  types  of  systems.  However,  one  must  keep  in  mind  that  if 
the  planning  curve  extends  over  many  thousands  of  hours,  the  planned  growth  rate  may 
not  be  sustainable  due  to  resource  constraints  besides  test  time  and  due  to  technological 
constraints.  Past  comparable  growth  rates  may  have  been  estimated  from  test  data  over 
one  test  phase  of  a  much  shorter  duration  and  the  system  may  also  have  been  relatively 
immature. 
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4.  DERIVED  RELIABILITY  GROWTH  PATTERNS. 


4.1  Assumptions. 

The  system  has  a  large  number  of  potential  failure  modes  with  initial  rates  of 
occurrence  .  The  modes  are  candidates  for  corrective  action  if  they  are  surfaced 

during  test.  All  failure  modes  independently  generate  failures  according  to  the 
exponential  distribution  and  the  system  fails  whenever  a  failure  mode  occurs.  It  is  also 
assumed  that  corrective  actions  do  not  create  new  failure  modes. 


4.2  Background  Information. 


The  first  step  in  obtaining  a  functional  form  for  the  expected  failure  intensity  as  a 
function  of  test  time  and  planning  parameters  that  is  based  on  non-empirical 
considerations  involves  the  relationship  between  the  expected  number  of  failure  modes 
surfaced  and  test  duration.  This  relationship  was  considered  by  Crow  (1982)  for  the  case 
where  test  duration  is  continuous.  In  this  paper  we  are  measuring  test  duration  in  a 
continuous  fashion.  Test  time  will  be  used  as  a  generic  measure  of  test  duration  for  this 
continuous  case.  The  relationship  is  easily  obtained  by  expressing  the  number  of 
surfaced  modes  by  test  time  t  as  a  sum  of  mode  indicator  functions.  In  particular,  let 
L(t)  denote  the  indicator  function  for  mode  i.  The  indicator  function  takes  on  the  value 
one  if  mode  i  occurs  by  t  and  equals  zero  otherwise.  The  number  of  modes  surfaced  by  t 
is  given  by, 

M(0  =  XC(0  (5) 

i= 1 

The  expected  value  of  M(t)  is  equal  to, 


k  k  k 

m = J  £(/,«) =X(l“eH,) = k  ~ 

/=!  i=l  i= 1 


(6) 


This  expected  value  function  implies  a  functional  form  for  the  expected  failure  intensity 
and  corresponding  MTBF  as  a  function  of  test  time  t  given  that  corrective  actions  have 
been  incorporated  to  all  the  failure  modes  surfaced  by  t.  One  component  of  the  expected 
failure  intensity  is  due  to  the  failure  modes  not  yet  surfaced  by  t.  This  component  is 
simply  given  by  the  derivative  of  jLi(t) .  Note, 


(7) 


In  (Ellner  et  al.,  2000)  it  is  shown  that  the  expression  in  (7)  is  the  expected  failure 
intensity  due  to  all  the  modes  not  surfaced  by  t.  To  show  this  observe  that  the  failure 
intensity  due  to  these  modes  can  be  expressed  as  the  random  variable  A  u(t)  where 


Auto  =  5>M.-(')) 

i=l 


The  expected  value  of  At/  (r)  is  given  by, 


(8) 


k  k 

E{Av(t))=  £  2,.  {t -£(/,.(;))}  =  J  V* 

1=1  /=! 


djLi(t) 

dt 


(9) 
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Before  considering  the  other  components  of  the  system  failure  intensity,  we  shall 
address  obtaining  parsimonious  approximations  to  the  expected  number  of  failure  modes 
surfaced  by  t  and  the  corresponding  failure  intensity  due  to  the  unsurfaced  failure  modes. 
The  exact  expressions  for  these  quantities  are  given  by  (6)  and  (9).  Note  that  these 
expressions  depend  on  k  +1  parameters,  namely  the  number  of  potential  failure  modes  k 
and  the  initial  failure  mode  rates  of  occurrence  At  for  /  =  1. . 


4.3  Parsimonious  Approximations. 


4.3.1  Expected  Number  of  Modes  and  its  Derivative. 


To  obtain  parsimonious  approximations  to  the  expected  number  of  modes  surfaced  by 
t  and  its  derivative,  we  shall  consider  an  optimization  problem  under  the  assumption  that 
all  corrective  actions  are  delayed  until  t.  Let  Nt  denote  the  number  of  failures  that  occur 

by  t  due  to  mode  i.  Then  A.  =  —  denotes  the  standard  Maximum  Likelihood  Estimate 

t 

(MLE)  of  Ai .  Consider  the  estimator  for  A,  given  by, 


Ai=6-Ai+(\-6)-avg(Ai)  (10) 

where  avg{Xi)  denotes  the  arithmetic  average  of  the  k  i,  and  0  e  [0.1]  is  chosen  to 


minimize  the  expected  sum  of  squared  errors  between  Ai  and  At ,  i.e.  E 


£(4-4)2 


.  The 


value  of  9  that  solves  this  optimization  problem  can  be  shown  to  be  9S  (Ellner  et  ah, 
2004)  where 


0S  = 


Var[A 


A 

k-t 


1~!  1  +  VarlAJ 


(11) 


k 

for  A=YjA,  , 

i= 1 


and 


Var[Aj  \  = 


Ew-*)2 

1=1 _ 

k 


The  estimate  of  At  given  by  (10)  with  9 


equal  to  0S  has  been  called  the  Stein  estimate  of  At  (Ellner  et  ah,  2004).  Note  this  is  a 
theoretical  estimate  in  the  sense  that  it  cannot  be  computed  from  the  data  since  it  involves 
the  unknown  values  of  k,  X,  and  Var[Af] .  The  quantity  Var[A ,.]  can  also  be  expressed  as 
follows: 


Var\ A,  ]  = 


(12) 


From  the  definition  of  A,  and  the  fact  that  Ni  equals  zero  for  a  failure  mode 
unobserved  by  t  we  have  that  the  Stein  assessment  for  the  failure  rate  contribution  of  a 
failure  mode  not  observed  by  t  is  given  by, 

(13) 

Thus  the  Stein  assessment  of  the  failure  intensity  due  to  all  the  failure  modes  not  surfaced 
by  t  equals, 
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!*=(*-»)•( !-*,)■£  =  1-y  -(1-^)-  y 


where  m  denotes  the  number  of  surfaced  modes  by  t  and  obs  denotes  the  index  set  for  the 
failure  modes  not  surfaced  by  t.  From  (11)  and  (12)  one  can  show, 


i -es  = 


•Ml 

U -t) 

i  t 

A2  \  A 
-  +  - 

k  I  k-t 


Replacing  1  -6S  in  (14)  by  the  final  expression  in  (15)  and  simplifying  yields, 


1  \  /  \  f  \T 

j  1  j  m  N 

k  J  1  k  )  \  t 


K 

n  I* 

1--  +  - 


Equation  (16)  gives  the  Stein  assessment  for  the  failure  intensity  due  to  all  the  failure 
modes  not  surfaced  by  t.  Note  that  in  can  be  regarded  as  an  estimate  of  the  expected 
number  of  modes  surfaced  by  t,  i.e.  ju(t) .  Additionally,  in  light  of  Equation  (9),  the  left 

hand  side  of  (16)  can  be  viewed  as  an  estimate  of  ^-1 .  From  (7),  it  follows  that  the 

dt 

derivative  of  /./(f)  at  t  =  0  equals  X.  Finally,  observe  that  —  in  (16)  is  the  maximum 
likelihood  estimate  of  the  initial  failure  rate  X  under  the  assumption  that  all  corrective 

actions  are  delayed  to  t.  Let  h(t)  denote  ^—1 .  Simulation  results  for  a  number  of  cases 

dt 

(where  k  and  the  2,.  are  known)  conducted  in  support  of  (Ellner  et  ah,  2004)  have 
indicated  that  the  Stein  assessment  given  in  (16)  yields  good  estimates  of  h(t)  when  all 
the  corrective  actions  are  delayed.  The  value  h(t )  that  is  being  estimated  does  not  depend 

on  the  corrective  action  process.  Only  the  estimate  of  X  given  by  —  depends  on  the 
assumption  that  all  corrective  actions  are  delayed  until  t.  Thus  the  right  hand  side  of  (16) 
with  m  and  —  replaced  by  good  approximations  of  //(/)  and  A  =  c^!l  respectively 


with  m  and  —  replaced  by  good  approximations  of  //(/)  and  A  =  c^!l  respectively 

t  dt  l=0 

should  yield  a  good  approximation  for  h(t )  irregardless  of  the  corrective  action  process 
for  the  cases  where  the  Stein  estimate  of  h(t)  given  by  (16)  are  accurate.  This  suggests 
that  we  choose  our  parsimonious  approximation  for //(f) ,  denoted  by  juk(t),  as  the  unique 

solution  to  the  differential  equation  obtained  from  (16)  by  replacing  in  by  juk(t) ,  y  by  X  , 


and  by  with  initial  conditions  juk( 0)  =  0  and  dfUk^ 


=  A .  For  the  case  where 


all  the  2,  are  equal,  one  can  show  that  //,  (/  )  =  //(/  )  for  all  t  >0  .  Thus,  in  what  follows  we 
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shall  only  consider  the  case  where  not  all  the  Xi  are  equal.  The  solution  to  the  resulting 
differential  equation,  for  this  case,  with  the  specified  initial  conditions  is, 


Mk(t)=k\ l-(l  +  A  -t)~Pk  _ 


(17) 


where 


A 


x_ 

k 


(18) 


and 


Pk  = 


X 


(19) 


The  solution  was  obtained  by  the  method  of  integrating  factors  (Boyce  et  ah,  1965).  The 
solution  can  be  verified  by  directly  substituting  (17)  and  its  derivative  into  the  differential 
equation  for  juk(t)  and  noting  that  juk(t)  satisfies  the  specified  initial  conditions.  Observe 
that  jLik  (t)  can  be  expressed  in  terms  of  t  and  three  constants,  namely  k,  /.  and  pk .  The 


corresponding  parsimonious  approximation  for  h(t)  is 


^  which  we  shall  denote  by 
dt  : 


Kit). 


It  is  interesting  to  note  that  juk{t)  given  by  (17)  is  the  same  expression  that  one 
can  obtain  for  the  expected  number  of  software  bugs  surfaced  in  execution  time  t  given 
by  the  doubly  stochastic  exponential  order  model  presented  in  (Miller,  1985)  for  the  case 
where  the  initial  bug  occurrence  rates  Xl,...,Xk  constitute  a  realization  of  a  random 
sample  of  size  k  from  a  gamma  random  variable.  The  density  function  of  this  random 
variable  is  given  by, 


/(*)  = 


r(«i  +i)Al“‘*l) 


0 


for  x  >  0 
otherwise 


In  this  density  function,  T  denotes  the  gamma  function,  pk  is  defined  by  ( 1 8),  and  a,  + 1 
equals  pk .  This  result  is  shown  in  (Ellner  et.  ah,  2000)  where  pk  (t)  denotes  the  expected 
number  of  surfaced  modes  by  time  t  that  will  be  mitigated  by  a  corrective  action. 


4.3.2  Expected  System  Failure  Intensity  and  MTBF. 

Next  we  shall  consider  the  expected  system  failure  intensity  after  t  test  hours  and 
a  corresponding  parsimonious  approximation,  given  that  corrective  actions  are 
implemented  to  all  the  surfaced  failure  modes.  We  shall  let  df  denote  the  fraction 
reduction  in  the  rate  of  occurrence  of  mode  i  due  to  the  corrective  action  (termed  a  fix). 
The  reduction  factor  is  termed  the  fix  effectiveness  factor  (FEF)  for  failure  mode  i.  Let 
A (t)  denote  the  failure  intensity  of  the  system  given  that  fixes  have  been  incorporated  to 
all  the  failure  modes  surfaced  by  t.  Then, 

A(0  =  £(t-<W0K  (20) 

i=l 

The  corresponding  expected  failure  intensity  is  p(t)  where, 
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(21) 


pin  =  E{Mt))  =  £  (1  -  d,E{ /,.  «)K  =  £  (1  -  4  K  +  £ 

*=i  i=i  i=i 

This  expression  for  the  expected  failure  intensity  was  presented  in  (Crow,  1982). 

For  reliability  growth  planning  purposes,  assessments  of  individual  failure  mode 
FEFs  will  not  be  available.  Thus,  in  place  of  pit)  we  shall  use  a  parsimonious 
approximation,  denoted  by  pk(t) ,  that  utilizes  an  average  fix  effectiveness  factor.  It 
follows  from  (9)  that  A- hit)  is  the  expected  failure  intensity  due  to  the  failure  modes 
surfaced  by  t  prior  to  mitigation.  Assume  these  modes  are  mitigated  with  an  average  FEF 
of  pd .  Then  the  expected  failure  intensity  due  to  the  surfaced  failure  modes  after 
mitigation  can  be  approximated  by  (l  -  pd)-(A.~  hit)) .  Thus  the  parsimonious 
approximation  for  pit)  will  be  defined  as  follows: 

pA0  =  (i-//J-(2-M0)+M0  (22) 

We  also  define  the  parsimonious  MTBF  approximation  of  MTBFit)  =  (pit))  for  reliability 
growth  planning  by  MTBFk  it)  =  (pk  it))  . 

For  planning,  it  can  be  useful  to  add  a  term,  XA  ,  to  the  expressions  for  p(t)  and  pk(t) 
given  by  (21)  and  (22)  respectively.  This  term  represents  the  failure  rate  due  to  all  the 
failure  modes  that  will  not  be  corrected,  even  if  surfaced  ...  referred  to  as  A-modes 
(Crow,  1982).  This  term  for  planning  purposes  would  be  given  by  the  quantity 
(l  -MS)-X .  However,  since  this  term  does  not  contribute  to  the  difference  between  p(t) 
and  pk(t)  we  shall  not  consider  it  further  in  this  section  or  Section  5. 

It  may  be  difficult  to  select  a  value  of  k  for  planning  purposes.  For  complex 
systems  or  subsystems  it  is  reasonable  to  use  the  limiting  forms  of  juk(t),  hk(t ),  and  pk(t) 
as  k  oo  .  Consider  the  limit  as  &  ->  oo  of  these  functions.  In  taking  the  limit  we  shall 
hold  X  fixed  and  assume  the  limit  of  pk  is  positive  as  k  increases,  say  /?„  e  (0,  oo) .  Under 
these  conditions  one  can  show  the  three  functions  converge  to  limiting  functions  which 
we  shall  denote  by  p„(t),  h„(t),  and  pm(t),  respectfully.  One  can  show, 

Fc(t)=  -j-  ln(l  +  PK-t)  (23) 

and 

h„(t)  =  ^-  =  x(l  +  p„.ty  (24) 

at 

Also,  px(t)  is  given  by  (22)  with  hk(t)  replaced  by  hm(t). 
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5.  SIMULATION. 


5.1  Simulation  Overview. 

We  wish  to  compare  the  parsimonious  approximations  to  realized  and  expected 
reliability  growth  patterns  with  respect  to  a  number  of  quantities.  To  do  so  we  shall 
generate  a  number  of  realized  reliability  growth  patterns  via  simulation  in  Mathematica. 
We  shall  consider  cases  where  the  failure  mode  initial  rates  of  occurrence  are  realizations 
of  a  specified  parent  population  for  several  choices  of  the  parent  distribution.  We  shall 
also  generate  reliability  growth  patterns  for  a  deterministically  specified  sequence  of 
failure  mode  initial  rates  of  occurrence  that  have  been  found  to  be  useful  in  representing 
initial  bug  rates  of  occurrence  in  software  programs  under  development  (Miller,  1985). 
The  simulation  consists  of  the  following  steps: 

(1)  Specify  inputs.  This  includes  items  such  as,  1.  test  duration,  2.  the  number  of  failure 
modes,  and  3.  the  sequence  or  parent  population  governing  the  initial  mode  failure 
rates. 

(2)  Produce  mode  initial  failure  rates.  Failure  rates  are  either  stochastically  generated,  or 
deterministically  calculated.  In  the  stochastic  case,  failure  rates  are  generated  by 
drawing  realizations  of  a  random  sample  from  a  specified  gamma,  Weibull, 
lognonnal  or  loglogistic  (Meeker  et  ah,  1998)  parent  population.  In  the  deterministic 
case,  failure  rates  are  calculated  in  accordance  with  a  specified  geometric  sequence. 

(3)  Generate  mode  failure  times.  The  mode  failure  times  are  generated  via  a  function  of 
randomly  generated  uniform  numbers,  and  the  mode  initial  failure  rates. 

(4)  Generate  mode  fix  effectiveness  factors.  The  FEFs  are  generated  by  drawing 
realizations  of  a  random  sample  from  a  beta  distribution  with  mean  0.80,  and 
coefficient  of  variation  0.10. 

(5)  Examine  quantities  and  plots  of  interest. 

5.2  Simulation  Results. 

Results  below  display  plots  of  the  expected  and  realized  number  of  surfaced 
failure  modes  for  stochastic  2,  generated  from  a  loglogistic  (Figure  1),  and  deterministic 
2,.  calculated  from  a  geometric  sequence  (Figure  3).  Also  shown  are  plots  of  the 
reciprocals  (i.e.  MTBFs)  of  the  expected  and  realized  system  failure  intensities  for 
loglogistic  2,.  (Figure  2),  and  geometric  2,.  (Figure  4).  The  geometric  initial  mode 
failure  rates  are  given  by 

2,.  =a-b'  (25) 

for  i  =  l,...,k  where  0 < a  and  0<Z><1.  All  the  displayed  quantities  have  been  averaged 
over  ten  replications  of  simulation  steps  (2)  through  (4)  above. 

The  intent  of  the  plots  is  to  see  whether  the  functional  form  of  the  parsimonious 
approximations  are  reasonably  compatible  with  respect  to  (1)  the  expected  number  of 
surfaced  failure  modes  as  a  function  of  test  time,  and  (2)  the  reciprocal  of  the  expected 
system  failure  intensity  as  a  function  of  test  time.  Corrective  actions  are  assumed  to  be 
implemented  to  all  the  failure  modes  surfaced  by  t  with  the  simulated  mode  fix 
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effectiveness  factors. 


The  value  of  pd 


in  Equation  (22)  is 


1  2  , 

set  equal  to  — 


to 


generate  the  parsimonious  approximations  to  the  exact  expected  failure  intensity  and 
corresponding  MTBF.  Additionally,  for  the  results  displayed  below,  A- =  1,500  and 
A  =  10  1 . 

The  value  of  the  scale  parameter  obtained  from  Equation  (18)  does  not  provide 
adequate  parsimonious  approximations  except  when  the  parent  population  is  gamma  or 
the  scale  parameter  is  sufficiently  small.  Thus  for  the  specified  k,  jud ,  and  A  ,  the  scale 
parameters  J3k  and  p.,  of  the  parsimonious  approximations  were  fitted  to  the  exact 
expected  number  of  surfaced  failure  modes  function  by  using  maximum  likelihood 
estimates.  These  estimates  were  obtained  from  the  simulated  mode  first  occurrence 
times.  This  was  accomplished  by  assuming  the  generated  initial  mode  failure  rates 
AI,...,Ak  represented  a  realization  of  a  random  sample  of  size  k  from  a  gamma 

distribution  with  scale  parameter  pk  and  mean  ^ .  This  procedure  provided  a  “best 

statistical  fit”  of  the  parsimonious  functional  approximations  for  p{t)  and  pit) ,  with 
respect  to  the  scale  parameter,  over  the  entire  planning  period  of  interest,  i.e.  10,000 
hours. 

The  parsimonious  approximations  for  pit)  and  pit)  based  on  the  limiting  forms 
for  pk(t)  and  pk(t)  as  k  increases  will  tend  to  be  too  large  for  values  of  t  when  pit)  is  too 
close  to  k.  We  have  observed  that  the  limiting  approximations  are  adequate  for  pit)  and 
pit )  over  the  range  of  t  for  which  p(t)<lyi  .  Thus  for  complex  systems,  or  subsystems, 
the  limiting  approximation  functional  forms  should  be  adequate  representations  of  pit) 
and  pit )  over  most  test  periods  of  interest. 

The  red  curves  in  the  figures  below  represent  the  exact  expected  number  of 
surfaced  modes  (Figures  1,  and  3),  or  the  reciprocal  of  the  exact  expected  system  failure 
intensity  (Figures  2  and  4).  The  dots  in  each  figure  represent  a  corresponding  stochastic 
realization.  The  green  curves  display  the  finite  k  approximations  while  the  blue  curves 
display  the  corresponding  limiting  approximations.  The  displayed  curves  and  stochastic 
realizations  are  averages  over  ten  replications  of  simulation  steps  (2)  through  (4). 
Similar  results  were  obtained  for  the  cases  where  the  Ai  were  generated  from  gamma, 
lognonnal,  and  Weibull  parent  populations. 

For  comparison  purposes,  the  MIF-HDBK-189  system  MTBF  based  on  Equation 
(4)  was  fitted  to  the  reciprocal  of  the  expected  system  failure  intensity  (the  red  curves). 
The  MIF-HDBK-189  curves  are  displayed  in  yellow  and  were  fitted  utilizing  all  the 
observed  simulated  cumulative  times  of  failure.  The  use  of  all  cumulative  failure  times 
requires  that  fixes  be  implemented  when  failure  modes  are  observed.  The  simulation  was 
carried  out  in  this  manner  to  allow  the  parameters  of  the  MIF-HDBK-189  curves  to  be 
statistically  fitted  via  the  maximum  likelihood  estimation  procedure  in  (Department  of 
Defense,  1981).  As  for  the  other  displayed  quantities,  the  averages  of  10  replicated  MIF- 
HDBK-189  MTBF  curves  are  shown. 


Page  14 


700 

600 

“  500 

1 

»  400 

.  300 

£  200 

100 

- 

y 

0  2000  4000  6000  8000  10000 

Time 

Figure  1.  Avg.  No.  Surfaced  Modes  (Loglogistic). 


Notice  the  high  degree  of  accuracy  displayed  in  Figure  1  for  the  finite  and  infinite 
k  approximations  despite  violating  the  gamma  assumption  used  to  statistically  fit  the 
parsimonious  approximations. 


Figure  2.  Reciprocal  of  the  Failure  Intensity  (Loglogistic). 

Figure  2  displays  a  high  degree  of  accuracy  for  the  statistically  fitted  PM2  MTBF 
approximations  despite  violating  the  MLE  assumption  that  the  initial  mode  failure  rates 
are  gamma  distributed.  In  addition,  the  PM2  approximations  of  the  MTBF  appear 
favorable  to  that  of  the  MIL-FIDBK-189  model. 

Figure  3  and  4  below  are  analogous  to  Figures  1  and  2,  respectively.  The  only 
difference  is  the  generation  procedure  associated  with  the  initial  mode  failure  rates 
utilized  in  the  analysis.  In  this  case,  failure  rates  are  deterministically  calculated  in 
accordance  with  a  geometric  sequence.  The  results  are  similar. 


Figure  3.  Avg.  No.  Surfaced  Modes  (Geometric). 
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Figure  4.  Reciprocal  of  the  Failure  Intensity  (Geometric). 
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6.  USING  PLANNING  PARAMETERS  TO  CONSTRUCT  THE  PARSIMONIOUS 
MTBF  GROWTH  CURVE. 


6.1  Methodology. 


6.1.1  Planning  Formulae  not  Using  Failure  Mode  Classification. 


In  the  previous  sections,  a  functional  form  for  the  planned  MTBF  growth  curve 
was  developed.  It  was  indicated  that  this  functional  fonn  was  compatible  with  a  number 
of  potential  growth  patterns.  In  Section  5  the  simulation  produced  failure  mode  first 
occurrence  times  from  a  set  of  initial  mode  failure  rates.  For  each  simulation  replication, 
the  parsimonious  MTBF  growth  pattern  was  derived  from  a  statistically  fitted 
parsimonious  expression  for  the  expected  number  of  failure  modes  function.  This  was 
accomplished  by  utilizing  the  mode  first  occurrence  times  to  obtain  a  MLE  of  scale 
parameter  [i  subject  to  the  initial  failure  intensity  X  held  fixed  to  a  specified  value  (e.g.  X  = 
0.10  in  Section  5).  In  practice,  the  initial  mode  rates  of  occurrence  will  not  be  available 
to  obtain  the  planning  curve  parameter  [>. 

In  this  section  we  shall  develop  formulas  for  /?  in  terms  of  the  planning  parameters 
T,  Mh  Mg,  and  average  FEF  pd  (and  k  for  the  finite  case).  We  shall  also  address  the 
question  of  how  well  the  parsimonious  MTBF  planning  curves  based  on  the  resulting 
values  of  [)  captures  several  potential  reliability  growth  patterns  that  depend  on  realized 
values  of  Xl,...,Xk  and  dl,...,dk. 

As  indicated  in  Section  4,  the  form  of  the  parsimonious  expected  system  failure 
intensity  is, 

PpL(t)=(1-Md\X'-h{t))+h(t)  (26) 


For  complex  systems, 


h(t) 


X 

1  +  j3  •  t 


(27) 


where  0 <t<T .  For  the  finite  k  case,  the  equation  for  pPL(t)  remains  the  same  with  h(t) 


replaced  by  hk(t)  = - — r  where  pk  =  ^A.  p  and  pk  denotes  the  planning  value  of /? 

(1  +  Pk  ■  t)  k 

for  finite  k. 

To  develop  formulas  for  [1  in  terms  of  planning  parameters,  let  0(t)  denote  the 
expected  fraction  of  X  attributed  to  the  failure  modes  surfaced  by  t.  Thus, 

0^)=X-h(A  =  l_hA)  (28) 

A  A 

This  yields, 

h(t)=x-(\-e(t))  (29) 

It  follows  that, 

Ppl (0=  (1  - ){X  ■  e{t)}+x  •  {1  - 8{t)\=\\  - pd  •  0{t)}- X  (30) 


Let  Mg  denote  the  goal  MTBF  at  t  =  T  and  Xa  =  MGl .  Then  we  set 

ka  =  PP,{T)  =  [\-Pd-e(T)\-X 


Thus, 


(31) 
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(32) 


0(T)  = 


Pa 


1-- 


M, 


M 


a  J 


For  finite  k  let  0{t)=  0k (t)  where  0k(t)=  =  pk  -t) 

A 

solution  to  the  equation  (32)  with  d(r)  =  0k (t)  where  p,  =  • 

Note  for  the  complex  system  case, 

e(t)=x-^  15 -1 


A  1  +  (3  •  t 


Therefore,  for  this  case 


Solving  for  P  yields, 


p-T 
l  +  P-T 


1_ 

VVd) 


X_M, 


V  MGJ 


p= 


1- 


Mr 


Md 

V  v 


\_M^ 


M 


oJJ 


In  the  above  f3k  is  the 


(33) 

(34) 


(35) 


6.1.2  Planning  Formulae  Using  Failure  Mode  Classifications. 


In  some  cases  the  set  of  failure  modes  can  be  split  into  two  categories  tenned  A- 
modes  and  B-modes  (Crow,  1982).  The  B-modes  are  failure  modes  that  will  be  mitigated 
if  surfaced  during  test.  The  A-modes  are  those  that  will  not  receive  a  corrective  action 
even  if  observed  during  test.  For  this  case,  the  parsimonious  expected  failure  intensity 
would  be, 

P PL  (0  =  2^}  +  (l  —  Pd  \^-B  —  ^fl(0)+^7s(0  (36) 

where  XA  is  the  failure  intensity  due  to  A-modes,  AB  is  the  initial  failure  intensity  due  to 
B-modes  (thus  A  =  Aa+Ab),  hB(t )  is  the  expected  failure  intensity  due  to  the  set  of  B- 
modes  not  surfaced  by  t,  and  ju'd  is  the  average  FEF  that  would  be  realized  for  the  B- 


modes  if  all  were  surfaced  during  test.  For  complex  systems,  hB(t)  is  given  by  - — . 

Using  an  argument  similar  to  the  one  in  Section  6.1.1  it  can  be  shown  that  for  this  case 
planning  fonnula  (35)  becomes, 


P  = 


M  r 


MS  •  n'd  - 


t-^ 

.  Mc  JJ 


(37) 


where  MS 


The  planning  parameter  MS  is  tenned  the  management  strategy.  This 


represents  the  fraction  of  A  that  is  due  to  the  initial  B-mode  failure  intensity.  For  the 


finite  k  case,  hB(t)  is  given  by 


An 


(1  +  ft-0P'k 


P'k+l 


where 


Pk  = 


k-Pk 


The  value  ph  solves 
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Equation  (32)  with  6{t)  replaced  by  6ti :B  (7)  =  l-(l  +  f3k -7)  ^Pk+^  and  pd  replaced  by 
MS  ■  p'd  . 


6.2  Comparisons  of  MTBF  Approximations  Using  Planning  Parameters. 


In  what  follows,  we  shall  not  use  failure  mode  categories.  Unlike  the  planned 
MTBF  growth  curve,  MTBFPL  (?)  =  pP[  (?) ,  the  average  MTBF  growth  path  generated  from 
the  simulation  replications  depends  on  the  particular  parent  population  of  the  2,  or 
deterministic  formula  used  to  generate  the  2,  ,  together  with  the  generated  mode  FEFs 
drawn  from  a  beta  distribution.  Thus,  this  average  MTBF  growth  path  over  ?e[o,r] 
depends  on  far  more  than  just  k,  T,  M  f ,  MG ,  and  pd .  Hence  one  cannot  expect  that  the 


planned  growth  path  from  M ,  to  MG ,  based  solely  on  the  planning  parameters,  will 
always  closely  match  the  averaged  reciprocals  of  the  exact  expected  system  failure 
intensity.  However,  as  indicated  in  the  preceding  sections,  the  functional  fonn  of  the 
parsimonious  MTBF  planning  curve  is  more  compatible  with  respect  to  the  realized 
MTBF  growth  pattern  than  the  MIL-HDBK-189  power  law  MTBF  growth  pattern. 
Additionally,  the  planning  parameters  are  easier  to  interpret  and  directly  influence  than 
those  utilized  in  the  MIL-HDBK-189  approach. 

In  a  number  of  instances  of  practical  interest  the  parsimonious  MTBF  model 
based  on  the  planning  parameters  closely  approximates  the  averaged  exact  MTBF  growth 
patterns.  To  consider  this,  we  shall  compare  the  parsimonious  MTBF  planning  curve  to 
the  reciprocals  of  the  realized  stochastic  system  failure  intensity  and  expected  system 
failure  intensity.  For  a  given  simulation  replication  we  shall  stochastically  generate  from 
a  given  parent  population  or  deterministically  calculate  2, , . . . ,  A ,  together  with 
corresponding  mode  FEFs  dx,...,dk.  The  FEFs  are  generated  on  each  replication  from  a 


beta  distribution  with  a  mean  of  0.80  and  coefficient  of  variation  of  0.10. 

To  calculate  the  planning  value  of  fd  on  each  simulation  replication,  we  shall  set 


k 

M7=2~’  where  2  =  ^2.  and  choose  pd 

i=l 


(one  could  alternately  choose  jud  to  be 


the  expected  value  of  the  beta  distribution).  The  value  of  Mc  is  set  equal  to  the 
reciprocal  of  the  realized  value  of  the  stochastic  system  failure  intensity  at  t  =  T.  Then 
Equation  (32)  with  the  appropriate  form  of  0(t)  is  used  to  obtain  the  planning  [j  for  the 
finite  k  and  complex  system  cases.  The  corresponding  finite  k  and  complex  system 
MTBF  planning  curves  for  the  replication  are  given  by  MTBFPL  (?)  =  pp[  (?)  where  pPL  (?)  is 
specified  in  Equation  (26). 

The  plots  below  in  Figures  5,  6,  8,  and  10  compare  the  average  of  ten  replicated 
MTBF  finite  and  infinite  k  planning  curves  (green  and  blue  curves,  respectively)  to  the 
corresponding  averages  of  the  reciprocals  of  the  following  failure  intensities:  (1) 
stochastic  realizations  of  the  system  failure  intensity  (black  dots);  (2)  the  expected  system 
failure  intensity  (red  curves)  and;  (3)  MIL-HDBK-189  planning  curve  failure  intensities 
(yellow  curves).  For  our  examples  the  test  period  is  7=  10,000  hours  and  k  =  1,500. 

One  problem  encountered  in  utilizing  planning  parameters  to  generate  the  MIL- 
HDBK-189  curves  is  that,  as  noted  in  Section  3,  these  curves  employ  an  average  MTBF 
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over  a  selected  initial  test  phase  (since  the  curves  interpolate  back  to  a  zero  MTBF  at  t  = 
0).  To  generate  a  MIL-HDBK-189  planning  curve  on  each  simulation  replication  using 
Mj ,  Mg  ,  and  T  we  have  used  the  initial  system  MTBF  metric  given  in  (Crow,  2004)  for 
power  law  growth  curves.  The  expression  for  M,  given  in  (Crow,  2004),  denoted  by 
MI  Ce  ,  is  given  by 


M 


I,CE 


(  n  \^Pm 

\Am) 


(38) 


where  XM  and  pM  are  the  MIL-HDBK-189  power  law  parameters  utilized  in  Equations 
(1)  and  (2).  The  rationale  for  the  M,  metric  given  by  (38)  is  not  discussed  in  (Crow, 
2004).  However,  it  may  have  been  motivated  by  the  following  fact:  for  a  non- 
homogeneous  Poisson  process  of  the  number  of  failures  experienced  by  test  time  t  with 
mean  value  function  2M  -tPu  ,  the  time  to  first  failure  is  Weibull  distributed.  Moreover, 
the  mean  of  this  Weibull  random  variable  is  equal  to  the  right-hand  side  of  Equation  (38). 

Averages  of  ten  replicated  power  law  MTBF  approximations  shown  in  yellow  in 
the  figures  below  were  obtained  using  the  M ,  metric  presented  in  (Crow,  2004)  as 
follows: 


1. 

2. 

3. 

4. 


k 

Set  M,  =2  1  where  2  =  ^2,.  and  21,...,2Jt  are  the  realized  failure  mode  initial 

i=l 

rates  of  occurrence  for  the  replication; 

On  each  simulation  replication  of  the  realized  stochastic  system  failure 
intensity,  set  Mc  equal  to  the  reciprocal  of  the  realized  failure  intensity  at  /  = 


T; 

Set  Mj  =  M ICE  and  M G  =  [lM  •  pM  ■  TPu  1 }  .  Solve  for  2M  and  pM  ,  and; 

Set  MTBF{t)={pM(t\l  for  0 <t<T  where  pM(t)=A.M  ■  PM  •  ( />M  1 .  The  values  for 
2m  and  pM  that  satisfy  the  two  equations  in  step  3  can  readily  be  shown  to 


satisfy  the  equations  M,  =  - 


fr 

v  Pm  )\ 

,  and  2m  = 


kPm-MgJ 


Pm 


j'^-Pm 

Pm 


First  we  consider  several  cases  where  the  2,.  are  considered  a  realization  of  a 
random  sample  of  size  k  from  a  specified  parent  population.  The  reference  (Miller,  1985) 
considers  a  class  of  software  reliability  models  where  the  initial  bug  rates  of  occurrence 
are  assumed  to  represent  such  realizations. 


6.2.1  Gamma  Parent  Population. 


Page  20 


40 

n  30 
H 

E 

4J 

u  20 
cd 

X 

w 

10 

******** 

_  _  « 

_ " 

0  2000  4000  6000  8000  10000 

Test  Time 

Figure  5.  Reciprocal  of  the  Failure  Intensity  (Gamma). 


In  Figure  5  above  the  parent  population  is  a  gamma  distribution.  This  distribution 
has  been  utilized  as  a  parent  population  for  initial  bug  rates  of  occurrence  in  the  software 
reliability  literature  (Fakhre-Zakeri  et.  ah,  1992).  For  this  case  the  averages  of  the  MTBF 
planning  curves  for  the  finite  and  infinite  cases  are  quite  close  to  the  averages  of  the 
reciprocals  of  the  realized  stochastic  and  expected  system  failure  intensities.  The  close 
approximation  is  a  consequence  of  the  relation  mentioned  in  Section  4  between  the 
solution  of  the  differential  equation  that  defines  the  parsimonious  function  juk(t)  for  the 
expected  number  of  failure  modes  surfaced  by  t  and  the  gamma  distribution. 

6.2.2  Lognormal  Parent  Population. 


Figure  6.  Reciprocal  of  the  Failure  Intensity  (Log  Normal). 


In  Figure  6  above  a  lognonnal  parent  population  is  utilized.  For  this  population, 
the  averaged  parsimonious  MTBF  approximations  based  on  the  finite  and  infinite 
planning  values  of  P  are  again  close  to  the  averages  of  the  reciprocals  of  the  realized 
stochastic  and  expected  system  failure  intensities.  However,  this  close  agreement  does 
not  always  occur,  as  illustrated  in  Figure  8.  To  consider  this  further,  we  shall  look  at  the 
portion  of  the  initial  failure  intensity  that  the  top  w  failure  modes  comprise  as  a  function 
of  w.  The  top  w  failure  modes  refer  to  a  set  of  failure  modes  of  size  w  whose  initial 
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failure  rates  are  at  least  as  large  as  the  initial  failure  rates  of  the  remaining  k-w  inodes. 
More  formally,  for  i  =  l,...,k  let  A(/j  denote  the  ordered  mode  initial  failure  rates  such  that 

Ay  >  •  •  •  >  .  Also  let  ti{w)=\  —  \ for  w  =  \,...,k  where  A  =  '^Al  . 


100  200  300  400  500  600  700 

Number  of  Modes 

Figure  7.  Top  W  Modes  (Log  Normal). 


The  graphs  of  the  average,  rj  (w),  of  the  ten  replicated  functions  77(11)  is  displayed 
in  Figure  7  for  the  case  where  on  each  replication,  A1,...,At  was  drawn  from  a  lognormal 
distribution  with  mean  •  This  lognormal  distribution  was  also  used  to  generate  the 

ten  replicated  graphs  of  the  realized  and  expected  system  failure  intensities  on  which 
Figure  6  is  based.  Note  in  Figure  7  that  77(100)  is  less  than  0.40. 
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Figure  8.  Reciprocal  of  the  Failure  Intensity  (Log  Normal). 


Figure  9  is  the  corresponding  graph  of  77  (w)  versus  w  for  the  MTBF  average 
curves  displayed  in  Figure  8.  The  ten  replicated  system  failure  intensity  curves  utilized 
in  Figure  8  are  based  on  the  same  ten  sets  of  initial  mode  failure  rates  used  to  construct 
the  graph  in  Figure  9.  As  for  Figures  6  and  7,  the  lognormal  parent  population  used  had  a 
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mean  equal  to  .  However  the  variance  was  larger  than  that  of  the  lognormal  utilized 
for  Figures  6  and  7.  This  resulted  in  77(100)  >  0.50  . 


In  general,  the  closer  the  graph  of  77  (w)  is  to  the  line  through  the  origin  with  slope 
equal  to  1  Ik,  the  closer  the  averaged  MTBF  planning  curve  will  be  to  the  averaged 
reciprocals  of  the  realized  stochastic  and  expected  system  failure  intensities.  This 
follows  from  the  fact  that  the  solution  to  the  differential  equation,  //,(/),  converges  to 
ju(t)  as  Ak  approach  a  common  value  20>0.  For  the  case  where  lt  =  A0  for 
7  =  1,..., k  on  each  replication,  the  graph  of  rj(w )  versus  w  for  w  =  l,...,k  lies  on  the  line 
through  the  origin  with  slope  equal  to  1/A:.  Note  the  graph  of  rj  ( vv)  in  Figure  7,  although 
not  close  to  the  line  through  the  origin  with  slope  1/A;,  is  closer  to  this  line  than  is  the 
graph  of  77  (w)  displayed  in  Figure  9.  This  is  consistent  with  the  planning  approximation 
in  Figure  6  being  more  accurate  than  in  Figure  8.  The  lognormal  distribution  for  the  2,. 
on  which  Figures  8  and  9  are  based  would  not  be  a  realistic  candidate  in  many  instances 
for  the  parent  population. 

6.2.3  Geometric  Initial  Mode  Failure  Rates. 

Finally,  we  consider  geometric  failure  rates,  a  case  where  the  initial  mode  failure 
rates  are  specified  deterministically.  For  this  case  recall  A.  =  a-b‘  for  i  =  l,...,k  where  a  > 
0  and  0  <  b  <  1.  According  to  Miller  (1985),  such  bug  initial  failure  rates  have  been 
estimated  during  replicated  -  run  software  debugging  experiments  (Nagel,  1984,  1982). 
For  such  a  deterministic  case,  only  the  set  of  associated  FEFs  dx,...,dk  are  regenerated 
from  the  beta  distribution  on  each  replication.  Note  that  one  can  show, 
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Thus, 


Also,  as  before  on  each  simulation  replication,  Aa  is  set  equal  to  the  realized  value  of  the 
stochastic  system  failure  intensity  at  t  =  T  and  Ma  =  Aj .  To  calculate  the  MTBF  planning 
curve  on  the  replication,  in  addition  to  T,  Mi,  and  Mg,  the  average  FEF  value  pid  must  be 


specified. 


We  set 


(an  alternate  possibility  would  be  to  set  jud  equal  to  the 


specified  mean  of  the  beta  distribution).  Figure  10  displays  the  averages  of  the  ten 
MTBF  planning  curves  for  the  finite  and  infinite  cases  over  the  replications  (green  and 
blue  curves,  respectively).  These  average  curves  are  compared  to  the  averaged 
reciprocals  of  the  realized  stochastic  system  failure  intensities  (black  dots)  and  the 
expected  system  failure  intensities  (red  curve).  Also  displayed,  is  the  average  of  the  ten 
MIL-HDBK-189  planning  curves  (yellow  curve).  These  are  based  on  T,  M G ,  and 


M  /  —  M I  CE  . 


Figure  10.  Reciprocal  of  the  Failure  Intensity  (Geometric). 


Figure  1 1  displays  the  graph  of  tj{w)  =  tj(w)  versus  w. 
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Figure  11.  Top  W  Modes  (Geometric). 


Even  though  the  top  100  modes  account  for  60%  of  the  initial  failure  intensity,  the  finite 
and  infinite  averages  of  the  10  replicated  PM2  planning  curves  are  reasonably  compatible 
with  the  averages  of  the  reciprocals  of  the  10  corresponding  realized  stochastic  system 
failure  intensities  and  expected  system  failure  intensity  graphs.  The  parameter  b  governs 
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the  percent  of  the  initial  failure  intensity  contributed  by  A{l),...,X{k).  One  can  show 


r/(w) 


\-b" 
1  -bk 


for  w  =  1, 


k  .  Also,  holding  X  and  b  fixed  we  obtain  iif  (w)  =  \[mri(w)  =  \-bw 

£->oo 


for  w  =  1,2, ... . 

As  for  the  lognormal  case,  if  77(100)  is  lowered  from  its  current  value,  then  this 
will  generally  lead  to  the  finite  and  infinite  average  of  the  PM2  MTBF  planning  curves 
being  closer  to  the  averages  of  the  reciprocals  of  the  realized  and  expected  system  failure 
intensities. 
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7.  RELIABILITY  GROWTH  POTENTIAL. 


7.1  Growth  Potential  in  Terms  of  Planning  Parameters. 


In  contrast  to  the  MIL-HDBK-189  planning  model,  the  PM2  planning 
approximation  recognizes  a  ceiling  for  the  MTBF.  To  obtain  this  limiting  value  we  note 
that  lim  h(t)  =  0  .  From  Equation  (26)  we  then  obtain 

t— >00 

lim  pPL{t)=(l- jud)-A  (41) 

00 

This  limiting  value  of  pPL(t)  is  termed  the  growth  potential  failure  intensity  and  is 
denoted  by  pap.  The  growth  potential  MTBF  is  defined  to  be  p(-,,{t)  and  is  denoted  by 


MCP.  Note, 


M  GP 


M, 

1  -  Fd 


(42) 


From  Equation  (30)  it  is  clear  that  pPL(t)  is  a  strictly  decreasing  function  of  t  whose  limit 
as  ?  ->•  oo  is  pGP .  Equivalently,  MTBFPL  (t)  =  p~[  (t)  is  a  strictly  increasing  function  of  t 
whose  limit  as  t  -» oo  is  MaP . 

The  above  comments  with  respect  to  pPL(t)  and  MTBFPL(t)  also  apply  to  the  case 
where  failure  modes  are  classified  into  A-modes  and  B-modes.  However,  for  this  case, 


M  GP  —  ■ 


M, 


i  -  (ms)- p'd 

where  p'd  now  denotes  the  average  FEF  with  respect  to  the  B-modes. 


(43) 


7.2  Planning  Parameter  /?  in  Terms  of  Growth  Potential. 


The  complex  system  planning  formula  for  /?,  given  by  Equation  (35),  can  be 
rewritten  in  terms  of  the  MTBF  growth  potential.  One  can  show, 


P  =  \- 

i  T 


_M, _ 

l  Ma 
M  Gp  J 


(44) 


This  formula  applies  whether  one  or  two  failure  mode  categories  are  utilized,  as  long  as 
the  appropriate  expression  for  MaP  is  applied.  For  a  logically  consistent  set  of  reliability 
growth  planning  parameters  one  must  have  M,  <  Ma  <  MaP .  Note  this  ensures  that  /?>(). 


7.3  Plausibility  Metrics  for  Planning  Parameters. 


Observe  that  Equation  (44)  shows  that  if  T  is  chosen  to  be  unrealistically  small  for 
growing  from  M ,  to  MG  then  the  resulting  value  of  f>  will  be  unduly  large.  This  would 
be  reflected  in  the  function, 


40 


P-t 
1  +  f5  •  t 


(45) 
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rising  towards  one  at  an  unrealistic  rate.  For  example,  a  large  [>  could  imply  that 
9(t0)=  0.80  for  an  initial  time  segment  [(U0]  for  which  past  experience  indicates  it  would 
not  be  feasible  to  surface  a  set  of  failure  modes  that  accounted  for  80%  of  the  initial 
failure  intensity.  An  unrealistically  large  /?,  and  corresponding  d(t)  function,  could  also 
arise  by  choosing  Mc  to  be  an  excessively  high  percentage  of  MCP .  This  discussion  also 
pertains  to  the  case  where  two  failure  mode  categories  are  utilized.  For  this  case  0(t )  is 
replaced  by  0B  (7) ,  or  the  left-hand  side  of  Equation  (45).  The  value  0B  (7)  denotes  the 


expected  fraction  of  AB  attributed  to  the  B-modes  surfaced  by  t.  The  scale  parameter  [> 
utilized  in  obtaining  0B  (7)  is  still  given  by  Equation  (44).  However,  to  obtain  0 B  (7) , 


Mgp  is  computed  via  Equation  (43). 

A  second  potentially  useful  metric  for  judging  whether  the  planning  parameters 
give  rise  to  a  reasonable  value  for  fi  is  the  implied  average  mode  failure  rate  for  the 
expected  set  of  surfaced  modes  over  a  selected  initial  reference  test  period  [o,t0]. 
Denoting  this  average  mode  failure  rate  by  pavg(t0)  we  have, 

<46> 

Classifying  failure  modes  into  A  and  B  modes  we  have, 


Pa 


(t  h B  (f)  ) 


(47) 


where  pavgJi  (t0 )  denotes  the  average  B-mode  failure  rate  for  the  set  of  B-modes  expected 
to  be  surfaced  during  [o,?0].  In  Equation  (47),  pB(t0)  is  the  expected  number  of  surfaced 
B-modes  over  [o,r0] .  For  complex  systems,  Equation  (46)  yields, 

(P-toY 


,(0=- 


(48) 


C  •(1+^-t0)-ln(1  +  A-to) 

where  P  is  given  by  Equation  (44)  with  M GP  expressed  by  Equation  (42).  Likewise  for 
complex  systems,  Equation  (47)  implies, 


p  (t  )  = _ ^,?°) _  (49) 

where  P  is  given  by  Equation  (44)  with  MGP  expressed  by  Equation  (44).  For  a  given  M , 
(and  MS  for  pmg,B(t0)),  Pmg{t0 )  and  pavgJ!(t(>)  can  be  expressed  in  terms  of  T  and  Mc  for 
M,  <Mg  <Map.  Denote  these  expressions  for  pavg(t0)  and  pavgJI(t,)  by  P,rJl\M G-t(])  and 
pmg  (r ,  M c ;  / o ) ,  respectively.  The  following  properties  of  these  functions  can  be  useful  in 
judging  whether  T  and  Mc  are  reasonable  planning  values: 

1.  pavg{T,Ma;t0)  and  pavgtB{T,MG',t0 )  are  continuous  positive  strictly  decreasing 
functions  of  T  and  strictly  increasing  functions  of  Ma  for  M ,  <MC  <  M ap ; 

2.  pmgir,M a;t0)  and  pavgtB{T,M a;t0)  approach  infinity  as  T  approaches  zero  or  Ma 


approaches  MCP  and; 

3-  Pavg{T’M g'Po)  and  pavgtB{T,MG-t0 )  approach  zero  as  T  approaches  infinity  or  MG 
approaches  M, . 
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A  third  potentially  useful  metric  to  judge  the  reasonableness  of  the  planning 
parameters  is  the  expected  number  of  unique  failure  modes  or  unique  B-modes  surfaced 
over  an  initial  test  interval  [(),  t(l  ]  they  imply.  Prior  experience  with  similar  developmental 
programs  or  initial  data  from  the  current  program  can  serve  as  benchmarks. 
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8.  GENERATING  A  PLANNED  RELIABILITY  GROWTH  PATH. 


Once  the  planning  parameters  are  chosen,  the  parsimonious  approximation  for  the 
expected  failure  intensity  can  be  used  to  generate  a  detailed  reliability  growth  planning 
curve.  For  example,  suppose  a  test  schedule  is  laid  out  that  gives  the  planned  number  of 
RAM  miles  accumulated  on  the  units  under  test  per  month.  Also  suppose  the  test 
schedule  specifies  blocks  of  calendar  time  for  implementing  corrective  actions.  Finally, 
for  planning  purposes  let  us  assume  that  in  order  for  a  failure  mode  to  be  addressed  in  an 
upcoming  corrective  action  period,  it  must  occur  four  months  prior  to  the  start  of  the 
period.  For  this  situation  the  MTBF  could  be  represented  by  a  constant  value  between 
the  ends  of  corrective  action  periods  and  between  the  start  of  testing  and  the  end  of  the 
first  scheduled  corrective  action  period  (CAP).  For  such  a  test  plan,  jumps  in  MTBF 
would  be  portrayed  at  the  conclusion  of  each  CAP.  The  increased  MTBF  after  the  jump 
is  given  by  MTBF(ti)={p(ti\x  where  ti  denotes  the  accumulated  test  time,  as  determined 
from  the  monthly  schedule,  by  the  calendar  date  that  occurs  four  months  prior  to  the  start 
of  the  CAP.  Since  MTBF (tt)  depends  on  a  large  number  of  parameters  it  would  be 
approximated  by  the  parsimonious  approximation  {/y(f )}  1  or  {pSh))  •  In  such  a 
manner  a  sequence  of  target  MTBF  steps  would  be  generated  that  grow  from  the  initial 
MTBF  value  to  a  goal  MTBF  value. 

Figure  12  below  depicts  a  detailed  reliability  growth  planning  curve  for  a 
complex  system  for  the  case  where  A  and  B  failure  mode  categories  are  utilized. 


Figure  12.  PM2  Reliability  Growth  Planning  Curve. 

The  blue  curve  in  Figure  12  represents  MTBF(t)={pPL(t)}  where  pPL(t)  is  given  by 
Equation  (36)  and  [j  is  obtained  from  the  planning  parameters  by  Equation  (37).  Note  the 
value  MTBF(t)  is  the  system  MTBF  one  expects  to  obtain  once  all  corrective  actions  to  B- 
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modes  surfaced  during  test  period  [o,t]  are  implemented.  The  MTBF  steps  are 
constructed  from  the  blue  curve,  the  schedule  of  corrective  action  periods  (CAPs),  and 
the  assumed  average  corrective  action  implementation  lag.  In  Figure  12  note  that  the 
goal  MTBF,  M  G ,  was  chosen  to  be  larger  than  M R  =  65  hours,  the  MTBF  to  be 
demonstrated  during  a  follow-on  Initial  Operational  Test  and  Evaluation  (IOT&E).  This 
test  is  an  operational  demonstration  test  of  the  system’s  suitability  for  fielding.  Such  a 
test  is  mandated  by  public  law  for  major  DoD  developmental  systems.  In  such  a 
demonstration  test  it  may  be  required  to  demonstrate  M R  with  a  measure  of  assurance. 
In  the  figure  we  have  utilized  as  our  measure  of  assurance  a  demonstration  of  M R  at  the 
80%  statistical  confidence  level.  To  have  a  reasonable  probability  of  achieving  such  a 
demonstration,  the  system  must  enter  the  IOT&E  with  an  MTBF  value  of  M+r  which  is 

greater  than  M R .  The  needed  value  of  M+r  can  be  determined  by  a  well-known  statistical 

procedure  from  the  IOT&E  test  length,  the  desired  confidence  level  of  the  statistical 
demonstration,  and  the  specified  probability  of  being  able  to  achieve  the  statistical 
demonstration.  After  determining  M+r  one  can  consider  what  the  goal  MTBF,  M G , 
should  be  at  the  conclusion  of  the  development  test.  The  value  of  Ma  should  be  the  goal 
MTBF  to  be  attained  just  prior  to  the  IOT&E  training  period  that  proceeds  the  IOT&E. 
The  goal  MTBF  associated  with  the  development  test  environment  must  be  chosen 
sufficiently  above  M+r  so  that  the  operational  test  environment  does  not  cause  the 

reliability  of  the  test  units  to  fall  below  M+r  during  the  IOT&E.  The  significant  drop  in 

MTBF  often  seen  could  be  attributable  to  operational  failure  modes  that  were  not 
encountered  during  the  developmental  test.  In  Figure  12  a  derating  factor  of  10%  was 

M  + 

used  to  obtain  M G  from  M+r  ,  i.e.,  in  the  figure  Ma  =  . 
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